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NONLINEAR DYNAMICS OF A DIFFUSING INTERFACE 


Walter M.B. Duval 

National Aeronautics and Space Administration 
Glenn Research Center 
Cleveland, Ohio 44135 

Excitation of two miscible viscous liquids inside a bounded enclosure in a microgravity 
environment has shown the evolution of quasi-stationary waves of various modes for a range of parameters. 
We examine computationally the nonlinear dynamics of the system as the interface breakup and bifurcates 
to resonance structures typified by the Rayleigh-Taylor instability mechanism. Results show’ that when the 
mean steady field is much smaller than the amplitude of the sinusoidal excitation, the system behaves 
linearly, and growth of quasi-stationary waves occurs through the Kelvin-Helmholtz instability mechanism. 
However, as the amplitude of excitation increases, nonlinearity occurs through subharmonic bifurcation 
prior to broadband chaos. 


Introduction 

We examine computationally the nonlinear 
dynamics of flow fields, which destabilizes an initially 
stable sharp interface between two miscible liquids inside 
a cavity. The twx> fluids are subjected to oscillatory 
vibration in a microgravity environment, and serve as a 
model to understand effects of g-jitter. The destability of 
the interface between the two miscible liquids caused by 
bifurcation of the flow field serves as a model problem to 
study transport phenomena in a wide range of 
microgravity applications such as crystal growth of 
acousto-optic optoelectronic materials, protein and 
solution crystal growth. We consider a parametric space 
of relevance to microgravity applications. 

The microgravity environment is particular 
useful for transport phenomena involving 

nonhomogeneous bulk fluids due to minimization of 
buoyancy induced flows. However, vibration effects, 
which cause g-jitter, can serve as a source to induce 
convective flows. In order to quantify effects of g-jitter 
we designed an experiment to study effects of prescribed 
sinusoidal forcing on the stability of two miscible liquids 1 . 
The experimental results show that for the frequency band 
of interest (0.1 - 10 Hz) and disturbance amplitude (0.2 - 
20 milli-g) the interface can indeed become unstable and 
generate quasi-stationary waves up to four modes. The 
instability of the interface was shown to cause by Kelvin- 
Helmholtz (K-H) instability 2 . The flow field destabilizes 
to an oscillating parallel shear flow and causes the 
interface to evolve like a vortex sheet. The parametric 
space for the experiment ranged from Stokes-Reynolds 
number (Res) of 0.002 to 0.5. In this computational study, 
we extend Res to 43. 


We address the dynamical state of the flow' field 
to gain insight on the predictability of the system. Of 
interest is the number of modes that can be generated for a 
quasi-stationary wave prior to its breakup. We examine 
the characteristic response of the system to various 
amplitude and frequency of excitation. We show that for 
Res O(.l) the flow field responds isochronously to the 
input excitation, thus generates a limit cycle. The number 
of modes generated at the interface is correlated with the 
magnitude of Res. With increasing Res, the system 
becomes nonlinear and responds with harmonics of the 
input frequency. A subharmonic bifurcation occurs, 
similar to a two torus bifurcation, prior to broadband 
chaos. In the neighborhood of the chaotic state, there is 
breakup of the interface. The mechanism of the breakup is 
due to Rayleigh-Taylor instability. 

In the following, we consider an infinite 
dimensional system as a dynamical model. We analyze the 
dynamical state of the system from its time history, 
pseudo and phase space trajectories, and power spectrum. 
The dynamical state of the system is correlated with the 
local bifurcation of the interface. 

Formulation 

The problem of interest is shown in Figure 1, two 
miscible liquids with distinct density meet initially at a 
planar interface. This realization has been shown to be 
possible in a microgravity environment. The two fluids are 
subjected to a body force parallel to its interface. The 
body force 

g(t) = ngo + mg 0 (cos co t +9) ( 1 ) 
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consists of the steady background DC component (ng 0 ), the left and right fluids respectively, with a value of 0.5 
the excitation AC component with amplitude (mgo), for the interface, 
circular frequency co , and phase angle 0 . 


The dynamical model of the problem consists of 
the Navier-Stokes equations coupled with species 
continuity. For a square cavity (H=L) with characteristic 
length, time and velocity scales selected as, 




1.0 

<05 

0 


0<jc<L/2 

x = L/2 


(9) 


L c =Vv/« , T=l/co, U c = (2) 

p CO 

The dimensionless incompressible Boussinesq 
equations in vorticity-stream function relationship can be 
expressed as, 


The concentration field equation is normalized with C* = 
(C-C B )/(C A -C B ). The parametric space (A) of the above 
set of equations 

A = A(Res, AR , Sc, 0 ) (10) 
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consists of the Stokes-Reynolds (Res) number, the 
amplitude ratio (AR), the, Schmidt number (Sc), and the 
phase angle 0 defined as, 
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The velocity components and vorticity are defined as, 

dy/ dw „ 8v du 

„=_^, v=- — , £=- —(6) 

ay ax ax ay 


The boundary conditions are no-slip along the wall 
boundaries (F) 
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For the parameters of interest, the Schmidt 
number is kept constant as well as the phase angle. Thus 
the flow field as well as the dynamical evolution of the 
interface can be expressed as a co-dimension two 
problem, 


V - V( x,y,t; Res, AR) (12) 


We used finite difference to solve the coupled set 
of equations (3-5). Since we are interested in the 
asymptotic dynamics of the system, we used an explicit 
scheme with a third order Adams-Bashforth technique for 
time discretization. The Flux Corrected Transport 
algorithm is used to resolve the sharp discontinuity of the 
interface 3 . A 90 x 90 grid size is used for spatial 
resolution and found to be adequate over the parametric 
range of interest. We mostly varied Res for our parametric 
studies, however, we found that the results are sensitive to 
the magnitude of AR. 


V • n — 0 on F (7) 

and the condition of impermeability normal {ft) to the 
boundary, 

VC • n = 0 on r (8) 

Since we are considering an initial value problem, the 
concentration field has prescribed values of 1 and 0 for 


Discussion and Computational Results 

The parametric space considered consists of 
variation of the amplitude of excitation (mgo), input 
frequency (to), while keeping the background DC 
acceleration (ng o =10’ 6 go, go=980cm/sec 2 ) and phase angle 
(0 =0) fixed. The geometric length scales (H=L=5cm) and 
properties of the fluid (Sc=1079) are fixed. The range of 
parameters considered is shown in Table 1 . 
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Table 1 - Range of Parameters, ngo=10‘ 6 go 


f 

(Hz) 

Res 

AR 

0.5 

0.1908 

2.5 xi O' 4 

0.5 

0.2386 

2.0 xlO -4 

0.5 

0.2863 

1. 67x10"* 

0.5 

0.3818 

1.25xl0 4 

0.5 

0.4772 

1.0x10"* 

0.5 

0.7158 

6.67x1 0' 5 

0.5 

0.95 

5.0x1 0' 5 

0.33 

1.78 

2.50X10' 5 

0.25 

2.69 

2.50x1 0' 5 

0.5 

2.85 

1.67xl0" 5 

0.5 

3.80 

1.25X10' 5 

0.5 

4.75 

1.0x10"”’ 

0.5 

5.70 

8.33xlO' 5 

0.5 

6.67 

7.14X10" 6 

0.5 

7.60 

6.25xl0" 6 

0.5 

8.55 

5.56x1 O' 6 

0.5 

9.5 

5.0xl0" 6 

0.5 

10.45 

4.55 xl O' 6 

: 



0.5 

42.75 

l.llxlO" 6 


The metric used to quantify the response of the interface 
to various excitation input is its wavelength (X) which is a 
function of Res and AR, 

X = URes, AR) (13) 

Even the co-dimension two bifurcation problem is a 
formidable computational task because the system loses 
and gains stability over a dense parametric set similar to 
the logistic map model 4 . Even though all the fine details of 
the bifurcation sequence are not computed, the range of 
parameters in Table 1 gives a broad glimpse on typical 
scenarios. Over most of the parametric space, we varied 
mgo, a narrow range of co was considered. 

Effect of Flow Field on Dynamics of Interface 

In order to correlate the dynamical equilibrium of 
the flow field in phase space with the interface response, 
we show in Figure 2 the local bifurcation of the quasi- 
stationary wave as a function of Res. The wavelength of 
the q-s wave decreases, or the number of modes increases, 
as Res increases. Breakup of the interface occurs in the 
neighborhood of Res=8.55 for which q-s waves no longer 
exist. The number of modes is dictated by the number of 
bends in the interface, an S shape has two modes. 
Transition from two to four mode q-s wave occurs for 


Res>0.29. Typical interface response of four mode q-s 
wave is shown for Res=0.72; the four mode response 
extends up to Res=1.78. The mode number at the interface 
is sensitive to the threshold of AR, for example at 
Res=0.34 and AR=5 x 10' 5 a four mode structure can be 
suppressed by increasing AR to 1.25 x 1CT 4 , only two 
modes result. There is a transition from 4 to 8 primary 
modes for Res=2.85; note 3 secondary modes occur with 
very small wavelength at the top and bottom boundaries. 
This bifurcation occurs in the neighborhood of Res=1.9. 

The mechanism which generates q-s waves for 
Res 0(2.85) is basically the Kelvin-Helmholtz (K-H) 
instability. A row of parallel vortices is generated along 
the interface. These vortices serve as rotors, which stretch 
and fold the interface. As Res increases, Res 0(4.75), 
there occur secondary waves riding on the primary modes 
on the interface. The K-H instability mechanism manifest 
itself perpendicular to the interface to form vortex rows. 
The interface no longer forms q-s waves, rather stretches 
and folds continuously to generate roll-up structures as 
typified by the Rayleigh-Taylor instability 5 . Only the 
initial state of the interface structure is shown, the flow 
field becomes quite complex with vortices occupying the 
entire cavity. The state of the flow r field is chaotic as will 
later be shown. 

The wavelength (X) of the interface, up to the 
breakup, decreases asymptotically as Res increases, 
shown in Figure 3. The wavelength is estimated from its 
average value at the center of the cavity. The dependence 
of X on Res follows a power law 

X c< C Res (14) 

C is a constant and n is less than 1; the decrease in 
wavelength is similar to a turbulent cascade. In contrast 
the maximum magnitude of the velocity , shown in 

Figure 4, increases with Res and also follows a power law 
in which n is greater than 1 , 

V^C Res * (15) 

The functional relationship of X and V max on Res illustrate 
the global trends. Much more details would be needed to 
show the bifurcation details as found for example in the 
logistic map. 

Dvamical State of the Flow Field 

We characterize the dynamical state of the system 
from the time history of the asymptotic dynamics of the 
flow field. The limit set of the time history is obtained 
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from both the pseudo and phase space trajectories. For a 
given frequency and amplitude of excitation, the response 
of the system is obtained from its power spectrum. The 
estimated power spectrum for the v component of the 
velocity field is calculated from, 




\l*V(t)e- i2Kf 'dt 


( 16 ) 


in which R is the finite duration of the time interval. The 
smooth spectral estimate P v ( / ) is obtained using a 
Hanning window for filtering, and a convolution relation 
as, 


Pv(f)=Wjf)*? v (f) (17) 


similar to the Rossler funnel 8 . The flow' field responds 
with harmonics (1.0, 1.5 Hz) of the fundamental input 
frequency of 0.5 Hz; the 1.5 Hz component is very' minute 
and is barely discernible. The gradual decrease in the 
amplitude of the u component of the velocity field gives 
rise to the dense limit cycle attractor. The combination of 
the signal of the u and v component is shown in the (u, v) 
phase-space trajectory, which shows the evolution of a 
nonchaotic attractor. It is known that strange nonchaotic 
atttractors exist for quasi-periodically forced systems 9 ’ 10,11 , 
however they are not smooth manifolds in phase space. 
Even though the attractors shown in Figure 6 are 
nonchaotic, they are not strange. This study gives insight 
into parametric regions for which an infinite-dimensional 
system may display these so-called strange nonchaotic 
attractors. 


W ,(/) is the spectral window for filtering. The power 
spectrum is estimated using the Cooley-Tukey algorithm, 
which employs the use of Fast Fourier Transforms 6 ' 7 . The 
dynamics of the flow field is evaluated from a 
representative point of the flow field. For each data set we 
kept 10,000 points. 

We now analyze the dynamical state of the flow 
field, which generates the q-s waves shown in Figure 2. 
For low Res=0. 1 9 which yields a two-mode q-s w r ave, the 
flow field oscillate isochronously, see Figure 5a, with the 
input frequency f=0.5 Hz. The power spectrum shows that 
the v component of the flow 7 field responds with a 
frequency of 0.5 Hz, the same scenario was found for the 
u component. The pseudo-phase space trajectory, found 
by lagging the velocity with a constant k, of the v 
component indicates that the dynamics approach a limit 
cycle attractor for the selected location (0.33, 0.56) in the 
flow field. This indicates that for the given parametric 
input to the system, its future behavior will always tend to 
oscillate in phase with the input frequency. 

When nonlinearity becomes important. Figure 5b 
Res=0.48, harmonics of the fundamental frequency 
appear. While the u component of the velocity field is 
periodic, as shown above, the v component of the velocity 
shows a decrease in amplitude with every oscillation. The 
combination of the signal is shown in the phase space 
trajectory, w'hich shows a nonchaotic attractor. This 
attractor is a smooth surface with the attributes of 
stretching and folding in phase-space. 

With increasing nonlinearity in the system. Figure 6 
Re$=0.72, the flow field responds with a modulated 
amplitude signal in the v component. Such modulation 
gives rise to nonchaotic attractors as evidence from the 
pseudo-phase space trajectory of v(t + k); this attractor is 


A rational subharmonic bifurcation is shown in 
Figure 7. When Res is increased to 1.78 for an input 
frequency of 0.33 Hz, both the u and v component of the 
velocity field show harmonic response of the fundamental 
frequency. A complex nonchaotic attractor illustrating 
stretching and folding arises in phase-space. A 
subharmonic bifurcation occurred in the neighborhood of 
Res=2.69, AR = 2.5 x 10 5 . In contrast to the above cases, 
both the u and v components show response at one-half 
the input frequency of 0.25 Hz as well as higher 
harmonics. The response frequency of the system is 
phase-locked to the input frequency due to the rational 
relationship of the frequency ratio. Note that a 
subharmonic attractor is generated in phase space with a 
contraction of volume, since the system is dissipative; this 
reduces the complexity of the attractor in comparison to 
the lower values of Res. The amplitude ratio AR plays a 
significant role in the bifurcation sequence, for the same 
Res value of 2.69 if AR increases to 1.7 x 10 5 then the 
subharmonic bifurcation dissapears, only harmonics of the 
fundamental frequency occurs. 

The dynamics of the system for which it 
transitions from q-s waves to tendril structures and 
breakup is shown in Figure 8. Corresponding to the tendril 
structure, Res=4.75, is an irrational frequency component 
of the fundamental input frequency shown in the power 
spectrum. Similar to Res=2.69, there is a subharmonic 
attractor with contraction of volume in phase space. This 
corresponds to the tw r o-torus bifurcation as proposed by 
Newhouse, Ruelle, and Takens 12 . The mechanism of K-H 
instability which causes secondary waves on the interface 
gives rise to the complexity of the flow 7 field which 
generates these waves. In the neighborhood of the breakup 
region, Res=7.6, there is a resulting broadband power 
spectrum distribution which shows that the system 
responds to all the frequency components within the 
interval of 1.0 Hz. Since the dynamics of the system 
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changes after breakup, we partitioned the signal of the v 
component into two intervals. Contraction of volume 
occurs in phase space for the first interval. The second 
interval from 160 to 400 secs shows the asymptotic 
dynamics of the system which yields a chaotic attractor. 
This attractor has the attributes of folding and stretching 
which typifies unpredictability. The results for Res in the 
interval of 7.6 to 42.75 show that there is a sequence of 
window for which the system gains stability to become 
periodic and loses stability to become chaotic similar to 
the logistic map. 

Summary 7 and Conclusions 

We considered the nonlinear dynamics of flow 
fields, which destabilizes an interface through the Kelvin- 
Helmholtz instability mechanism to cause quasi-stationary 
waves. The wavelength decreases asymptotically as Res 
increases up to the breakup region; the decrease in 
wavelength is similar to a turbulent cascade). We show 
that when the system becomes nonlinear there is a 
transition from two to four mode q-s waves which gives 
rise to a host of nonchaotic attractors due to harmonics of 
the fundamental frequency in the response signal of the 
velocity component. Bifurcation of the system from q-s 
waves to tendril structures introduces a rational 
subharmonic bifurcation. Chaos occurs through an 
irrational subharmonic bifurcation, which corresponds to a 
two-torus bifurcation; the interface breaks up and the 
Raleigh-Taylor type of instability is manifested. Though 
Res was not varied continuously to capture all the 
bifurcation details, the subset of the parametric space 
considered lead insight into the dynamical behavior of the 
system. 
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Concentration field Stream function field Vector field 

t = 1 60 s, mg 0 = 1 5 mg, AR = 6.67 x 1 0-5, Re s = 0.72, V max = 0.049 cm/sec 
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Figure 2. — Bifurcation of quasi-stationary waves showing increasing density of modes up to breakup 
(Re s = 8.55), f = 0.5 Hz, ng 0 = 1 0^go- 
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Figure 5.— Limit cycle (Re s = 0.19) and nonchaotic attractor (Re s = 0.48) due to aperiodic component of velocity field. 
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Figure 6. — Nonchaotic attractor due to periodic and amplitude modulated components of velocity field. 
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